Hardcore Magnons in the 5 = 1/2 Heisenberg Model on the Square Lattice 
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We propose a versatile approach to treat commonly arising constraints. It is illustrated for inter- 
acting magnons of the Heisenberg antiferromagnet on a square lattice. For systems of L x L sites 
a non-perturbative continuous unitary transformation (CUT) is used to derive an effective Hamil- 
tonian conserving the number of magnons. They are bosonic particles with a hardcore constraint 
which is captured by a local, repulsive interaction U. The limits U — ► oo and L — *■ oo are achieved 
by extrapolation. The residual spin gap Ai is smaller than 0.01J reflecting the gapless nature of 
the magnons. The one-magnon dispersion displays all known characteristics. 
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Strongly correlated systems in low dimensions dis- 
play a large variety of fascinating properties like high-T c 
superconductivity quantum antiferromagnetism, quan- 
tum ferromagnetism or charge ordering. Theoretically, 
these phenomena are investigated in simplified models 
like the Hubbard model, the Heisenberg model and the 
t-J model. These models are derived from ab-initio calcu- 
lations, from more elaborate models, for instance three- 
band Hubbard models, or they are motivated by phe- 
nomenological considerations. 

In the t-J model transitions to double occupancies are 
forbidden. This represents a hardcore constraint. An 
analogous constraint arises if the magnons in an ordered 
Heisenberg model of spin S are treated as bosons. At 
most 2S magnons may be present at each site. A reliable 
treatment of such constraints is still a great challenge 
to condensed matter theory. Straightforward approaches 
like mean-field approximations fail because of the large 
(infinite) energy scale of the hardcore repulsion. 

In the present work we describe an innovative approach 
to treat hardcore constraints and other strong interac- 
tions. This approach is exemplified by the antiferromag- 
netic Heisenberg model on the square lattice. The ground 
state of this model displays long-range Neel order with 
a sublattice magnetization m « 0.3 < 1/2 reduced by 
quantum fluctuations. Therefore, a description in terms 
of magnons, propagating spin flips, as elementary excita- 
tions is a natural starting point which has inspired many 
studiesiSiSiiiSiSiLSiS. These treatments are often based 
on expansions in the inverse spin 1/S or in the inverse 
dimension l/d. So these approaches can be problematic 
for low spin in two dimensions. 

Magnons are interacting bosons subject to strong lo- 
cal constraints like the hardcore property. The hardcore 
property has been dealt with by choosing Hamiltonians 
with appropriate projections which ensure that no un- 
physical states are reached in the course of the Hamilto- 
nian dynamicsifi^iSiii. Recent approaches replace the 
magnon by other particles which do not exist freely or 
might exist only at higher energiesi^*i^. The purpose of 



the present work is to show that a direct treatment of 
the local constraint as an additional interaction is pos- 
sible and leads to very good results. This success opens 
the route to treat also more complicated models like the 
t-J model for doped Mott insulators in an analogous way. 
The Heisenberg model on the square lattice reads 

H=J 1 £s i S j (1) 

(id) 

where J > denotes the antiferromagnetic exchange con- 
stant, Si is the vector operator of a spin S = 1/2 on site 
i, and stands for the summation over nearest neigh- 
bors. We take the Neel state with J,-spins on the A-lattice 
and f-spins on the B-lattice as reference state. This ref- 
erence state will be mapped to the ground state by a 
continuous unitary transformation (CUT). The elemen- 
tary excitations are local spin flips which flip the |-spin 
up on the A-lattice or the "["-spm up on the B-lattice. 
Let (a) be the corresponding usual bosonic creation 
(annihilation) operators. These magnon excitations are 
hardcore particles, i.e., there can be at most one magnon 
per site. We impose this constraint by a local repulsive 
magnon-magnon interaction^ 

H u = UY,44H*i- (2) 

i 

The hardcore constraint is recovered for U — > oo. 
In terms of the magnons, the Hamiltonian reads 

H init = Yl h a k + B i { a W- q + h - c ) 

q 

+ Vk,k>, q al +q al_ q a k a k , (3) 

k,k' ,q 

up to a constant. All energies are given in units of 
J. The one-particle couplings arc uj q = 2, B q = 
\ (cos^a) + cos(g r J/ )), and the two-particle couplings are 
V ktk >, q = U/N - l/(2A0[cos(fc; + q x - k x ) + cos{k' y + 
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q y - k y ) + cos(efe) + cos(q y )] = U/(2N) + f(k, k' , q). 
The Hamiltonian remains invariant under translations by 
one lattice constant a (set to unity) in spite of the Neel 
state as reference state because a spin rotation by 180° 
has been performed on one sublatticei£. The two-particle 
couplings are manifestly symmetric under the exchange 
k <-> k' and q <-* k' + q — k. Numerically, we deal with a 
finite sample of linear size L with TV = L x L sites and 
periodic boundary conditions. 

A CUTi2ii£ is used to derive an effective Hamiltonian 
H eS which conserves the number of magnonsa****. This 
is done by solving the flow equation 



d l H(l) = [r,(t),H(l)} 



(4) 



where I is a continuous auxiliary variable and rj(£) the 
infinitesimal anti-Hermitian generator. The initial con- 
dition is H(£ = 0) = H mit . In the present work, the 
CUT is realized in a self-similar, renormalizing fashion 
in momentum space. The commutators for all expres- 
sions are computed using the standard bosonic algebra 
in momentum space. At maximum quartic terms in the 
bosonic operators are kept. Higher terms involving six 
bosons are neglected after normal-ordering with respect 
to the bosonic vacuum. So the Hamiltonian remains in 
the self-similar form of a Hamiltonian of pairwise inter- 
acting bosons 



H {i) = J2[4 e)a l a q + B i e) ( a W- q +^ 

q 

+ E V ki, q a l'+q a i-q a k a k' 



(5) 
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k x ,k 2 ,k 3 



\ , (at a] a] a, ,, ,, + h.c. 

l,k 2 ,k 3 l fcj k 2 k 3 k 1 +k 2 +k 3 



The anti-Hermitian generator rj(£) is chosen to be 

= 5>f ( t L 9 -h.c.' 



(6) 
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The effective Hamiltonian H := H(i = oo) is charac- 
terized by uig 



eff ._ 



and by V^ q 



V k,k',q > tiie 



other terms have to vanish. 

Comparing the coefficients of the same terms on the 
left and on the right hand side of Eq. 21 yields the (high 
dimensional) set of differential equations to be solved, 
which is done numerically. We treat systems of up to 
14 x 14 sites. The transformation can safely be carried 
out as long as the non-diagonal part of the Hamiltonian 
is decreasing monotonically 22 . As a measure of the non- 
diagonal part we define the residual off-diagonality 



ROD 2 = ^ q ( Bq 1 + ^k t ,k 2 ,k 3 (Afci,fc 2 ,fc 3 ) 



(7) 



Unfortunately, the numerical treatment of the flow 
equations reveals a significant increase of the ROD: no 
H can be obtained. In looking for the reason for this 
failure wc must keep in mind that the ROD decreases as 
long as it is dominated by terms in which the change of 
the number of elementary excitations is correlated to the 
change of the energy as measured by the diagonal part of 
the HamiltonianSii^&Si. In our system this means that 
an increase in the number of magnons has to imply an 
increase of the diagonal energy. If states with an incre- 
mented number of excitations are lower in energy the 
CUT breaks down. So a negative gap implies the failure 
of the CUT, see e.g. Ref. |23 This happens also in the 
system at hand. There is a bound state of two magnons 
with such a high binding energy that its excitation energy 
is negative. 

Why does this happen? The spin rotation symme- 
try requires that there has to be a tightly bound state. 
Dealing with finite clusters there is no true symmetry 
breaking so that the ground state is a singlet and the 
elementary excitation is a gapped triplet with three de- 
generate states S z £ {— 1,0,1}. The description in terms 
of magnons, i.e. spin flips up or down, canonically pro- 
vides two of them, namely at S z € {—1, 1}. The triplet 
with S z = has to be found in the sector with at least 
two magnons, one spin flip up and one down. The en- 
ergy of this state is equal to the energy of an elementary 
magnon, i.e. the S z = triplet must be a bound state of 
two magnons. This argument unambiguously shows that 
a magnon description of quantum antiferromagnets in or 
close to the paramagnetic phase implies a very strong 
interaction between these magnons. 

The fact that we neglected higher interactions spoils 
the delicate balance of magnon motion and magnon- 
magnon interaction. The numerics revealed that the 
attractive magnon-magnon interaction is overestimated 
leading to a too strongly bound S z = state. To remedy 
this problem we reduce the interaction by hand (see be- 
low) so that the spin symmetry, namely the degeneracy 
between the one-magnon states at Ai 

the two-magnon bound state at A2 

stored, i.e. Ai = A2. Here w 2 ^ 8 is the energy of the 
two-magnon bound state at total momentum K which 
can be determined by standard numerics once H has 
been found. For K = (0, 0), we solve the secular equation 
H cB E q l) - A 2 J2 q A q \q)\- q) where \q) stands 

for a single magnon at momentum q. The S = 1 state 
with S z = displays the symmetry A q+ ^,n) — —A q . 

The attractive interaction is deminished by reducing 
initial nearest-neighbor attraction by the reduction factor 
< A < 1 



^=(0,0) and 

2mae 

: = w if=(o,o)> 1S re " 



V k , v , q - Vtu- = U/(2N) + Xf(k,k',q) 



(8) 



The case A = corresponds to the total omission of the 
attractive magnon-magnon interaction keeping only the 
repulsive hardcore interaction. The actual initial Hamil- 
tonian H lmt is recovered for A = 1. The value of A is 
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fixed to A c where A2 = Ai holds. The above analysis 
of the failure of the direct CUT turns out to be valid. 
As long as A ^ A c the ROD decreases monotonically for 
large £ and the flow equations converge. 

The effective Hamiltonians H eS (L,U) were obtained 
for systems up to L = 14 and U = 500. The flow equa- 
tions were integrated until the ROD fell below 10 -4 . The 
large value of the interaction U necessitates a careful nu- 
merical treatment of the flow equations for small values 
of the flow I < 1/U leading to a slowing down of the 
program. Thus parallelization is mandatory to be able 
to deal with large systems. 

Fig- HI depicts a generic example how the bound state 
energy A2 depends on the reduction factor A. At A c ~ 
0.82592 a clear intersection is discernible; the correspond- 
ing gap value is Ai(Z = 8) = A 2 (L = 8) = 0.25008. The 
fact that A2 < 2Ai (dashed-dotted curve in Fig.^) shows 
that a true binding phenomenon is observed for appre- 
ciable interaction A ?t 0.76. For small values of A the 
repulsive interaction parametrized by U dominates and 
the system shows a large gap of the order J. 
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FIG. 1: Example of the one-magnon gap Ai (circles) and the 
two-magnon bound state energy A2 (squares) as function of 
the factor A in Eq. |H| The solid lines are splines; the inset 
zooms at the intersection Ai = A2. 

The same procedure has been performed for many val- 
ues of L and U. All reduction factors A c range from 
0.81 to 0.83. Extrapolating them in 1/U and in 1/L 
leads to X C (L = 00, U — 00) = 0.825. The deviation of 
X C (L = 00, U = 00) from unity, i.e. rs 17%, is a first es- 
timate for the size of the truncation error. In Fig. [21 the 
corresponding gaps Ai(L,/7) are shown. Linear extrap- 
olation in 1/U yields the values at U — 00. We found 
that a subsequent extrapolation in 1/L does not work. 
But an extrapolation in 1/ \fL works obviously very well. 
Except for L = 4, all data points lie nicely on straight 
lines. This makes us confident to proceed on the basis of 
the l/VL extrapolation although we do not know of an a 
priori reason for this unusual scaling. We presume that 
it results from using the Neel state as starting point for 
the treatment of finite clusters. 



It is reassuring that the residual gap A(L = 00, U = 
00) is as small as 0.006J. We take this fact as evidence 
that we could restore the spin symmetry by the procedure 
of reducing the initial attractive interaction in Eq. (jSJ). 
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FIG. 2: The one-magnon gap Ai(L,U) as function of l/VL. 
The value at U = 00 stem from linear extrapolation in 1/U. 
Solid lines are linear inter-/extrapolations. 
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FIG. 3: Upper panel: quantum correction factor Z c of the 
spin wave velocity is extrapolated in 1/U and 1/L; for com- 
parisons, see main text. Lower panel: inter-/extrapolation of 
the dispersion at two points of high symmetry. The points 
represent results extrapolated linearly in 1/U to U = 00. 

Next, we study the spin wave velocity c. It is charac- 
terized by the correction factor Z c which quantifies the 
renormalization relative to the result of spin wave theory 
Since finite lattices are studied the de- 
termination of a group velocity at vanishing momentum 
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is not possible. Thus we exploit the well-known shape of 
the dispersioniii and derive the velocity from 

where we use implicitly that the spin gap Ai vanishes. 
The upper panel in Fig. [3] depicts the results and their 
extrapolations. Since the bosonic dispersion aj° ff appears 

generically as square root y/ A\ + f(q) with a smooth 
function f(q) with /((0,0)) = it is consistent to ex- 
trapolate finite values w° ff in 1/L when the vanishing 

gap Ai scales like 1/y/L. We obtain Z S (L = oo,U = 
oo) = 1.14 which is about 3% away from the best val- 
ues Z c = 1.17947 (third order spin wave theory^) and 
Z c = 1.178(1) (high order series expansion 9 ). In view 
of the error estimate of 17% on the basis of 1 — A c the 
deviation in the low energy part of the spectrum is small 
and the result shows that the CUT approach is a valid 
way to tackle hardcore constraints. 

At high energies it is known that the difference be- 
tween ^(5/2 W2) ano - w (cf 7r) ^ s difficult to find. The third 
order spin wave results hardly display any dispersion be- 
tween q — (7r/2,7r/2) and q = (0,71")^. But series ex- 
pansion and quantum Monte Carlo calculations clearly 
show an appreciable difference between the dispersion 
at q = (71-/2,71-/2) and q = (0, 7r)L2i2. One finds at 
q — (tt/2, it/2) a saddle point with n , z \ ~ 2.17 and 

at q — (0, 7r) the maximum with uf^^ w 2.39. In the 
lower panel of Fig. |3| the corresponding results in our 
approach are displayed. The appropriate extrapolations 
finally yield ^^/ 27r /2) ~ 1-84 and s» 2.23. Hence, 

the former value deviates by about 15% and the latter 
one by about 7%. These deviations do not surprise in 
view of the first estimate of 17% on the basis of 1 — A c . 



We conclude that the omission of terms beyond the two- 
magnon interaction implies a truncation error of about 
17%. This is an encouraging result because it is surely 
possible to include at least the dominant parts of the 
higher magnon terms for further improvement. 

Further support for the approach chosen comes from 
the fact that w?5/2 n /2) 1S significantly lower than ufo^) 
which agrees qualitatively with the results by series ex- 
pansion and quantum Monte Carlo. The quantitative 
difference, however, is overestimated by about a factor 2 
which is attributed to the truncation of higher magnon 
terms. But the results show that the approach proposed 
captures the essential physics also at higher physics in 
contrast to diagrammatic spin wave theory. 

In summary, we propose non-perturbative continuous 
unitary transformations to treat the interactions aris- 
ing from constraints. This is successfully illustrated for 
magnons in the Heisenberg quantum antiferromagnet on 
the square lattice. All qualitative aspects of the disper- 
sion are retrieved when the symmetry between elemen- 
tary magnons and bound pairs of magnons is restored. 
This finding demonstrates that magnons in an antiferro- 
magnet constitute a strongly correlated system. Quanti- 
tatively, the results agree with previous findings within 
17%. Presently, the CUT approach is less accurate than 
the highly developed series expansions'^. But its accu- 
racy can certainly be enhanced and it is a versatile tool 
which can be used for the constraints in a multitude of 
models, including also doped systems, because it is for- 
mulated in the standard form of second quantization. 
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